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Barrier  islands  that  overlie  a  compressible  substrate,  such  as  islands  in  deltaic  environments  or  those  that 
overlay  mud  or  peat  deposits,  load  and  consolidate  the  underlying  subsurface.  Through  time,  the  elevation 
and  aerial  extent  of  these  islands  are  reduced,  making  them  more  susceptible  to  future  inundation  and 
overwash.  Sand  washed  over  the  island  and  onto  back-barrier  marsh  or  into  the  bay  or  estuary  begins  the 
consolidation  process  on  a  previously  non-loaded  substrate,  with  time-dependent  consolidation  as  a 
function  of  the  magnitude  of  the  load,  duration  of  load,  and  characteristics  of  the  substrate.  The  result  is  an 
increase  in  the  overwash,  migration,  breaching,  and  segmentation  of  these  islands. 

This  research  developed  a  two-dimensional  (cross-shore)  numerical  model  for  evolution  of  a  sandy  barrier 
island  that  spans  durations  of  years  to  decades  as  a  function  of  erosion,  runup,  overwash,  migration,  and 
time-dependent  consolidation  of  the  underlying  substrate  as  a  function  of  loading  by  the  island.  The  model 
was  tested  with  field  data  and  then  applied  to  evaluate  the  effects  of  a  compressible  substrate  on  long-term 
barrier  island  evolution.  Results  illustrate  that  barrier  islands  overlying  a  compressible  substrate  are  more 
likely  to  have  reduced  dune  elevation  due  to  consolidation,  incur  overall  volumetric  adjustment  of  the  profile 
to  fill  in  compressed  regions  outside  the  immediate  footprint  of  the  island,  and  experience  increased 
overwash  and  migration  when  the  dune  reaches  a  critical  elevation  with  respect  to  the  prevalent  storm 
conditions. 

Published  by  Elsevier  B.V. 


1.  Introduction 

Barrier  islands  form  a  dynamic  coastal  boundary  for  bays,  estuaries, 
and  mainland  shores.  They  buffer  these  more  fragile  environments 
from  coastal  winds,  waves,  and  storm  surges  (Stone  and  McBride, 
1998),  provide  habitat  for  static  and  migrant  populations  (Moore  et  al., 
1990),  and  foster  quiescent,  reduced  salinity  habitat  for  evolution  of 
juvenile  species  (Courrat  et  al.,  2009).  Estuaries,  particularly  those  on 
deltaic  coasts,  represent  the  most  productive  ecosystems  in  the  world 
yet  they  are  the  most  threatened  by  anthropogenic  activities  (Edgar 
et  al.,  2000),  an  inability  to  expand  with  relative  sea  level  rise  (so- 
called  “coastal  squeeze,”  French,  2006),  and  disintegration  of  protec¬ 
tive  barrier  islands  (e.g.,  McBride  et  al.,  1995;  McBride  and  Byrnes, 
1997;  Penland  et  al.,  2005).  Approximately  12%  of  the  world's  open- 
ocean  coast  is  fronted  by  barrier  islands,  and  28%  of  these  islands  occur 
in  deltaic  systems  (Pilkey  and  Fraser,  2003).  The  benefits  and 
functioning  of  barrier  islands,  especially  those  in  deltaic  settings, 
are  threatened  by  reduced  sources  of  sand,  relative  sea  level  rise, 
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and  anthropogenic  activities.  Intervention  to  restore  barrier  islands 
through  placement  of  beach-quality  sand  from  an  external  source  has 
been  conducted  since  the  1930s  (Farley,  1923;  Marine  Board,  1995) 
and  continues  to  be  considered  in  increasingly  large-scale,  regional 
applications  (van  Heerden  and  DeRouen,  1997). 

Barriers  evolve  in  form  and  migrate  in  response  to  coastal 
processes,  sediment  availability,  and  geologic  setting  on  time  scales 
ranging  from  hours  to  decades  to  centuries.  On  time  scales  of  hours  to 
days,  cross-shore  processes  during  storms  can  erode  the  foreshore  and 
overwash  the  island  and  deposit  sediment  on  the  back-barrier  or  into 
the  bay  (e.g.,  Leatherman,  1979;  Kahn  and  Roberts,  1982;  Dingier  and 
Reiss,  1990;  Doughty  et  al.,  2006).  Over  seasons,  years,  and  decades,  a 
gradient  in  longshore  sand  transport,  change  in  regional  sediment 
supply,  and  adjacent  inlet  processes  can  erode,  accrete,  and  migrate  a 
barrier  island  (e.g.,  Penland  et  al.,  2005;  Morton,  2008).  On  geologic 
time  scales  ranging  from  decades  to  centuries,  processes  in  the  vertical 
dimension  such  as  eustatic  sea  level  change,  regional  down-warping 
or  uplift,  and  consolidation  of  sediment  may  contribute  to  the  long¬ 
term  evolution  of  coastal  morphology  (e.g.,  Storms  et  al.,  2002;  Stolper 
et  al.,  2005;  FitzGerald  et  al.,  2008;  Moore  et  al.,  in  press).  For  barrier 
islands  overlying  poorly-consolidated  sediment,  such  as  deltaic,  bay, 
estuarine,  and  peat  deposits,  consolidation  of  the  underlying  substrate 
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Public  reporting  burden  for  the  collection  of  information  is  estimated  to  average  1  hour  per  response,  including  the  time  for  reviewing  instructions,  searching  existing  data  sources,  gathering  and 
maintaining  the  data  needed,  and  completing  and  reviewing  the  collection  of  information.  Send  comments  regarding  this  burden  estimate  or  any  other  aspect  of  this  collection  of  information, 
including  suggestions  for  reducing  this  burden,  to  Washington  Headquarters  Services,  Directorate  for  Information  Operations  and  Reports,  1215  Jefferson  Davis  Highway,  Suite  1204,  Arlington 
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due  to  the  weight  of  the  island  can  accelerate  long-term  morphologic 
response  (Guber  and  Slingerland,  1981;  Gayes,  1983;  Dean,  1997; 
Bourman  et  al.,  2000;  Rosati,  2009). 

Deltaic,  bay,  estuarine,  and  peat  deposits  compress,  or  consolidate 
as  a  function  of  the  load  that  is  applied,  duration  of  loading,  and 
characteristics  of  the  substrate  itself.  River  deltas  experience  consol¬ 
idation  wherever  the  river  deposits  organics  and  fine  sediments,  such 
as  silt  and  clay.  Deltaic  systems  that  experience  accelerated  subsidence 
include  the  Mississippi  River,  U.S.A.  (Coleman  et  al.,  1998);  Rhine- 
Meuse  River,  The  Netherlands  (Berendsen,  1998);  Ebro  River,  Spain 
(Sanchez-Arcilla  et  al.,  1998);  Nile  River,  Egypt  (Stanley  and  Warne, 
1998),  the  Ganges-Brahmaputra  Rivers,  Bangladesh,  India  (Allison, 
1998);  and  the  Yangtze  River,  China  (Xiqing,  1998).  Compaction  of  the 
subsurface  can  also  occur  with  fine  sediment  and  peat  deposited  in 
estuaries  and  bays  (e.g.,  Bloom,  1964;  Kaye  and  Barghoorn,  1964; 
Cahoon  et  al.,  1995;  Long  et  al.,  2006;  Meckel  et  al.,  2007).  For  barrier 
islands  overlying  a  soft  substrate,  whether  the  subsurface  is  riverine, 
estuarine,  or  organic  in  origin,  the  weight  of  the  island  compresses  the 
subsurface  resulting  in  a  time-dependent  reduction  of  the  island 
elevation.  The  net  effect  is  an  increase  in  the  likelihood  for  overwash  of 
the  island  and  subsequent  migration.  New  washover  deposits  begin  to 
consolidate  previously  non-loaded  sediment,  thus  perpetuating  the 
morphologic  change  process. 

In  this  paper,  we  develop  a  two-dimensional  (2D)  model  for  barrier 
island  erosion,  overwash,  and  washover,  including  time-dependent 
consolidation  of  the  subsurface  as  a  function  of  local  loading  by  the 
barrier  island  (Fig.  1).  The  model  is  applied  to  test  the  hypothesis 
that  cross-shore  barrier  island  migration  and  volumetric  losses  are 
modified  by  the  characteristics  of  a  compressible  subsurface,  as  a 
function  of  the  weight  and  duration  of  loading.  First,  we  review  the 
state-of-modeling  for  barrier  island  morphologic  change  over  a 
compressible  substrate  such  as  deltaic,  bay,  or  estuarine  deposits. 
Next,  we  present  theory  of  model  development,  validate  the  con¬ 
solidation  routine  with  long-term  data  from  the  Mississippi  River 
Delta,  and  illustrate  application  through  comparison  with  data  from 
two  barrier  islands  in  Virginia,  the  only  cross-shore  barrier  profile  and 
commensurate  sediment  core  data  presently  available.  We  apply  the 
model  to  evaluate  the  hypothesis  and  demonstrate  how  a  compress¬ 
ible  substrate  modifies  barrier  island  morphologic  evolution. 

2.  Behavior  and  modeling  of  barrier  islands  that  overlie 
compressible  sediment 

Barrier  islands  overlying  compressible  sediment  or  peat  deposits 
have  been  observed  to  rapidly  erode,  rollover,  breach,  breakup,  and 


possibly  become  submerged  (e.g.,  Penland  and  Boyd,  1981;  Leather- 
man  et  al.,  1982;  McBride  et  al.,  1995).  Penland  and  Boyd  (1981) 
defined  three  stages  of  deltaic  barrier  island  formation  based  on 
evolution  of  barrier  islands  associated  with  the  Mississippi  River 
Delta.  After  a  mature  active  delta  was  abandoned  by  the  river,  Stage  1 
began  with  an  erosional  headland  that  fed  flanking  barrier  islands 
with  sand  that  had  been  reworked  from  the  mixed  deltaic  deposit. 
Over  time  (millennia),  subsidence  and  wave-induced  erosion  deplet¬ 
ed  the  source  of  deltaic  sediment.  Stage  2  consisted  of  a  transgressive 
(retreating)  sandy  barrier  island  arc.  Finally,  Stage  3  occurred  when 
erosion  and  subsidence  reduce  the  barrier  island  to  a  subaqueous 
inner  shelf  shoal.  Until  human  intervention  began  in  the  early  1900s 
(through  levee  construction  and  river  diversion),  this  cycle  repeated 
as  the  river  occupied  new  locations  or  former  deltas  and  provided  a 
new  source  of  sediment. 

Because  of  this  cycle  of  delta  formation  and  abandonment,  the 
Louisiana  barrier  islands  are  comprised  of  a  relatively  thin  layer  of  fine 
sand  that  was  reworked  from  the  abandoned  delta.  The  islands  overlie 
a  thick  deltaic  sequence  of  clay  and  silt  that  was  deposited  during  the 
mid-to-late  Holocene  by  the  river,  and  eventually  transgressed  over 
back-barrier  estuarine  deposits  (Coleman  et  al.,  1998).  Penland  et  al. 
(2005)  documented  long-term  (greater  than  100  years)  and  short¬ 
term  (less  than  30  years)  shoreline  change  in  Louisiana  as  —6.1  and 
—  9.4m/yr,  respectively.  The  rapid  erosion  of  Louisiana's  coast  is 
attributed  to  the  predominance  of  muddy  sediment,  thickness  of  peaty 
marsh  soils,  rapid  rates  of  subsidence,  and  frequency  of  hurricanes 
(Kuecher,  1994;  Penland  et  al.,  2005).  Without  a  source  of  littoral  sand 
in  the  regional  coastal  system,  and  with  rapid  subsidence,  barrier 
islands  in  Louisiana  have  ultimately  drowned  (e.g.,  Ship  Shoal,  Penland 
and  Boyd,  1981).  In  the  region  of  the  abandoned  LaFourche  Delta, 
Kuecher  (1994)  correlated  thicker  deltaic  sediment  with  the  highest 
rates  of  land  loss  as  compared  to  thinner  deposits.  Similarly,  Penland 
and  Ramsey  (1990)  found  that  local  rates  of  relative  sea  level  rise  were 
related  to  the  thickness  of  Holocene  sediment  for  the  Mississippi  River 
Delta  and  Chenier  plains. 

Examples  of  barrier  islands  with  similar  cross-sections  occur  along 
the  Delaware-Maryland-Virginia  (“Delmarva”)  Atlantic  coast.  Some 
of  the  barrier  islands  in  this  region,  for  example  in  southwestern 
Delaware  Bay,  and  Wallops,  Assawoman,  and  Metompkin  Islands, 
Virginia,  are  comparable  to  the  deltaic  barriers  in  Louisiana  in  that  they 
are  composed  of  a  thin  veneer  of  sand  over  a  compressible  substrate. 
However,  in  this  region  the  origin  of  the  substrate  is  lagoonal  mud  and 
marsh  deposits,  over  which  the  islands  have  overwashed  and 
transgressed  through  time  (Kraft  et  al.,  1979;  Leatherman  et  al., 
1982;  Oertel  and  Kraft,  1994).  The  Delmarva  lagoonal  and  mud 
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Fig.  1.  Increase  in  overwash  and  migration  of  sandy  barrier  island  with  consolidation  of  underlying  compressible  substrate. 
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deposits  are  of  a  similar  thickness  (order  of  15  m  or  less,  Halsey,  1978) 
as  the  deltaic  sediment  sequences  in  the  vicinity  of  modern  Louisiana 
barrier  islands  (approximately  10  to  30  m;  Kulp  et  al.,  2002,  their 
Fig.  7). 

Conceptual,  analytical,  and  numerical  studies  incorporating  sand 
and  mud  fractions  in  long-term  barrier  island  evolution  and  migration 
are  being  developed  and  improved  (e.g.,  Cowell  et  al.,  1995;  van 
Maren,  2005;  Stolper  et  al.,  2005;  Alfageme  and  Canizares,  2005; 
Campbell,  2005;  Campbell  et  al.,  2005a, b;  de  Sonneville,  2006; 
Campbell  et  al.,  2006,  2007;  Thomson  et  al.,  2007).  Many  studies 
have  discussed  compaction  of  peat  deposits,  either  due  to  autocom¬ 
paction  (compaction  due  to  self- weight)  or  compaction  due  to  sub¬ 
sequent  loading  by  estuarine  sediments,  in  lagoons,  estuaries,  and 
deltaic  settings  (e.g.,  Bloom,  1964;  Kaye  and  Barghoorn,  1964;  Cahoon 
et  al.,  1995;  Long  et  al.,  2006;  Meckel  et  al.,  2007).  Peat  is  much  more 
compressible  than  other  types  of  substrate  such  as  sand,  clay,  silt,  or 
mud.  However,  similar  to  sediment  substrate,  the  magnitude  of  peat 
consolidation  is  related  to  the  thickness  of  the  deposit  (Kuecher,  1994; 
Meckel  et  al.,  2007).  For  example,  Bloom  (1964)  measured  13  to  44% 
compaction  of  a  sedge-peat  deposit  due  to  loading  by  a  10  m  deposit 
of  estuarine  mud  in  Clinton,  Connecticut.  Long  et  al.  (2006)  deter¬ 
mined  that  rapid  compaction  was  a  primary  mechanism  driving 
coastal  change  for  a  coastal  marsh  in  southeast  England,  United 
Kingdom,  where  the  peat  surface  compacted  at  least  3  m  due  to 
loading  by  4  m  of  intertidal  mudflat  and  tidal  channel  sediments. 
Meckel  et  al.  (2007)  developed  a  compaction  model  for  deltaic  settings 
and  concluded  with  a  statement  pertinent  to  this  study:  “high  density, 
permeable  sediments  such  as  sand,  at  the  surface  (typically  considered 
relatively  stable)  can  be  associated  with  high  compaction  rates, 
especially  if  they  overlie  thick  peat  deposits.” 

Guber  and  Slingerland  (1981;  see  also  Leatherman,  1987)  and 
Gayes  (1983)  considered  the  compression  of  back-barrier  and 
estuarine  sediment  through  placement  of  dredged  material  and  by 
washover  sand.  Data  from  two  dredged  sediment  disposal  sites  placed 
on  the  back  marsh  of  Assateague  Island,  Maryland,  indicated  a  linear 
relationship  between  the  effective  pressure  (overburden)  and  subsi¬ 
dence  of  the  marsh  surface,  with  the  older  site  having  greater 
subsidence  due  to  the  longer  loading  time.  Lateral  plastic  flow  and 
diapirism  (extrusion  of  sediment  from  the  substrate  such  as 
“mudlumps”  of  the  Mississippi  River  Delta  system,  Morgan,  1951) 
are  possible  with  loading  by  barrier  island  sands  and  tidal  deltas,  for 
cases  in  which  the  pore  water  pressure  is  the  same  order  of  magnitude 
as  cohesion  in  the  substrate  (both  in  units  of  kPa;  Guber  and 
Slingerland,  1981).  Compaction  and  lateral  flow  of  sediments  in  the 
vicinity  of  barrier  islands  include  the  following  potential  conse¬ 
quences.  (1)  Washover  removes  sand  from  the  foreshore  or  dune 
system  and,  with  settlement  into  a  compressible  marsh  or  bay 
subsurface,  will  induce  more  losses  from  the  barrier  island  system 
rather  than  increasing  elevation.  Thus,  the  island  will  be  susceptible  to 
additional  overwash  in  the  same  region  and  this  additional  overwash 
would  increase  migration  and  breakup  of  the  island,  suggesting  that 
compaction  may  also  be  a  factor  in  migration.  (2)  Barrier  island  cross- 
section  and  resulting  geomorphology  may  be  influenced  by  the 
subsurface  characteristics,  especially  if  the  subsurface  were  non- 
homogeneous.  Variable  retreat  rates  for  barrier  island  systems  might 
be  related  to  subsurface  characteristics  (e.g.,  void  ratio,  permeability, 
yield  criteria).  Observations  of  washover  fans  on  Assateague  Island, 
Virginia  and  Maryland  showed  lower  elevations  at  the  distal  ends  of 
the  fan  as  evidenced  by  ponding  as  compared  to  the  adjacent  marsh 
surface  which  was  at  the  mean  water  level.  These  fans  also  indicated 
that  lateral  flow  of  the  subsurface  sediments  might  have  occurred,  as 
evidenced  by  arcuate  ridges  paralleling  the  distal  portions  of  the  fans. 
(3)  The  potential  for  settlement  must  be  known  for  characterization  of 
the  sediment  budget  of  the  barrier  island  system  and  washover  fans, 
because  sequestration  of  sand  into  the  substrate  would  represent  a  net 
loss  to  the  budget. 


Roberts  et  al.  (1994)  linked  consolidation  of  Holocene  sediments  in 
the  Mississippi  River  deltaic  system  to  the  thickness  of  the  deposits.  In 
turn,  the  thickness  of  these  deposits  is  correlated  with  the  location  of 
previous  fluvial  entrenchment  by  the  Mississippi  River  system.  Thus, 
with  knowledge  of  the  sediment  type  and  former  locations  of  fluvial 
entrenchment,  Kuecher's  (1994)  and  Roberts  et  al.'s  (1994)  work  lend 
information  with  which  to  estimate  potential  future  compaction  of  the 
substrate  as  a  function  of  the  magnitude  of  the  loading. 

Bourman  et  al.  (2000)  discussed  rapid  geomorphologic  changes 
that  have  been  observed  at  the  River  Murray  Estuary,  Australia  due  to 
eolian,  riverine,  tectonic,  tidal,  and  wave  processes,  as  well  as  changes 
in  eustatic  sea  levels  over  the  past  125,000  years  and  recent  human 
activity.  Of  pertinence  to  this  review  is  the  observation  that,  as  barrier 
islands  fronting  the  River  Murray  Estuary  have  migrated  landward 
over  the  past  3000  years,  they  have  differentially  loaded  plastic  mud 
in  the  lagoon  resulting  in  an  increase  in  height  up  to  10  m  above 
present  sea  level.  The  authors  discussed  that  differential  loading  of 
these  lagoonal  sediments  was  sufficient  to  explain  their  elevation,  but 
that  seismic  events  may  have  also  played  a  role. 

This  review  of  experiences  demonstrates  that  loading  of  the 
substrate  through  barrier  island  migration  can  influence  the  resulting 
morphology  change.  In  particular,  subsurface  characteristics  can 
contribute  to  future  morphology  change  of  a  barrier  island  system 
through  three  mechanisms:  the  cross-section  of  the  island  as  a 
function  of  subsurface  characteristics  and  loading;  the  feedback 
between  consolidation,  barrier  island  elevation,  and  subsequent 
washover;  and  the  apparent  non-conservation  of  sediment  due  to 
compression  of  the  substrate.  None  of  the  conceptual,  analytical,  or 
numerical  models  that  were  reviewed  is  capable  of  calculating  these 
processes.  The  links  between  the  loading  by  the  islands  on  a 
compressible  substrate,  the  magnitude  of  time-dependent  consolida¬ 
tion  as  a  function  of  the  load,  and  the  subsequent  morphologic  change 
and  evolution  of  a  barrier  island  has  not  been  previously  quantified. 

3.  Model  development 

The  two-dimensional  Migration,  Consolidation,  and  Overwash  (2D 
MCO)  model  was  developed  to  investigate  the  role  of  consolidation  on 
migration  of  sandy  barrier  islands  that  overlie  compressible  sedi¬ 
ments.  The  paradigm  in  model  development  is  that  storm  waves  and 
elevated  water  levels  provide  the  primary  forcing  for  island  erosion, 
washover,  and  cross-shore  migration  during  an  annual  storm  season 
(e.g.,  see  Ritchie  and  Penland,  1988;  Dingier  and  Reiss,  1990,  1995; 
Rosati  and  Stone,  2009).  During  the  remainder  of  the  year,  sand 
transport  processes  are  assumed  to  be  relatively  minor  in  the  cross¬ 
shore  morphologic  evolution  and  migration  of  the  island.  Although 
longshore  transport  processes  during  this  time  could  translate  the 
island  laterally,  it  is  assumed  that  the  cross-shore  profile  is  not 
significantly  modified.  Time-dependent  consolidation  of  the  underly¬ 
ing  substrate  is  calculated  as  a  function  of  the  loading  by  the  island  as  it 
evolves  over  an  annual  storm  season.  The  island  evolves  over  years  to 
decades  as  a  function  of  storm  surge,  wave  height  and  period, 
consolidation  of  the  substrate,  and  the  rate  of  eustatic  sea  level  change. 

Initial  conditions  are  defined  by  a  sandy  barrier  island  with  a  given 
cross-shore  profile  that  overlies  a  sediment  substrate  of  specified 
characteristics.  The  2D  MCO  model  calculates  erosion,  runup  over¬ 
wash,  or  inundation  overwash  depending  on  the  storm  conditions  and 
relative  elevation  of  the  barrier  island.  If  washover  of  the  island 
occurs,  the  barrier  migrates  into  the  bay  and  consolidation  occurs  due 
to  the  existing  and  any  new  loading  (if  migration  occurred  onto 
partially-consolidated  sediments)  (Fig.  2). 

3.1.  Wave  transformation 

The  2D  MCO  model  transforms  storm  waves  from  deep  to  shallow 
water  using  a  time  series  of  wave  height,  period,  and  direction;  or,  the 
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Fig.  2.  Flowchart  for  2D  MCO  model. 


model  can  randomly  generate  wave  and  surge  conditions  from  user- 
specified  averages.  Deep-water  waves  are  transformed  from  offshore 
measurements  to  breaking  conditions  using  linear  wave  theory,  in  which 
time-dependent  measurements  of  deep-water  wave  height,  H0(t),  are 
related  to  wave  height  at  breaking,  Hb(t),  by  (Dean  and  Dalrymple, 
1984,  p.  115) 

Hb(f)  =  QV5/H02(t)C0(t)CQS60(t) 


as  beach  profile  morphology  using  the  Convolution  Storm  Erosion 
Method  (Kriebel  and  Dean,  1993), 

£(,)  =  T{’“TTsexp(“lr)“rrj?|cos2"' + 

(4) 

in  which  the  maximum  potential  recession  through  erosion,  ETO,  is 
given  by  a  form  of  the  Bruun  Rule  that  incorporates  wave  setup, 


The  depth-limited  breaking  criterion  is  K=Hb(t)/db(t)  =  0.78,  in 
which  db(t)  is  the  depth  at  breaking,  and  the  deep-water  wave  height 
and  direction  are  given  by  H0(t)  and  0o(t),  respectively.  The  deep¬ 
water  wave  speed  is  given  by  C0(t), 


Q(0 


Lott) 

T(t) 


(2) 


where  T(t)  is  the  wave  period  and  L0(t)  is  the  deep-water  wave  length 
equal  to 


Lott) 


gTitf 

2tt 


(3) 


Eoo  =  (S(t)  +  0.068Hb(t)) 


Wb(t) 

B  +  db(t) 


(5) 


In  Eq.  (5),  the  storm  surge  is  S(t),  Wb(t)  is  the  width  of  the  surf 
zone  out  to  db,  and  B  is  the  berm  elevation.  In  Eq.  (4),  /3t  represents 
the  percentage  of  time  the  storm  is  effective  in  eroding  the  beach  and 
(j=Tt/TD,  where  TD  is  the  total  storm  duration.  The  percentage  of 
effective  erosion  time,  /3t,  is  calculated  by 


(V 


2  TfL 


(6) 


in  which  the  duration  of  effective  erosion,  Ts,  is  given  by, 


3.2.  Erosion 


nb  (t)1J7  db(t)  tan 

gO.5^3  (l+  B  +  db(t )  ) 


Barrier  island  erosion  and  deposition  offshore  occur  if  the  storm 
surge  plus  wave  runup  do  not  exceed  the  barrier  island  elevation 
(Fig.  3a). 

The  time-dependent  recession  of  the  berm  through  erosion,  E(t),  is 
calculated  based  on  the  duration  and  magnitude  of  the  storm  as  well 


in  which  g  is  the  acceleration  due  to  gravity,  and  j3(t)  is  the  beach 
slope.  The  width  of  the  surf  zone  is  calculated  as 
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E(t)  =  Rate  of  erosion 

d(t)  =  Depth 

qDR  (t)  =  Rate  of  runup  overwash 

S  (t)  =  Surge  elevation 

qDi  (t)  =  Rate  of  inundation  overwash 

R  (t)  =  Runup  elevation 

zR  (t)  =  Surge  plus  runup  relative  to  dune 

bh  (t)  =  Dune  elevation 

qsw(t)=Rate  of  swash  transport 

Fig.  3.  Terminology  for  erosion  and  overwash  calculations. 


in  which  A  is  the  equilibrium  beach  profile  parameter  and  can  be 
related  to  median  grain  size  or  sand  fall  speed.  The  equilibrium  beach 
profile  concept  was  first  developed  by  Bruun  in  1954  (Bruun,  1962; 
Komar,  1998)  based  on  beaches  in  Monterey  Bay,  California,  and  the 
exponential  value  confirmed  by  Dean  (1977)  in  analysis  of  more  than 
500  beach  profiles  along  the  U.S.  Atlantic  and  Gulf  coasts  and  has  since 
been  applied  to  beaches  around  the  world  (e.g.,  Dean  et  al.,  1993; 
Larson  et  al.,  1999;  Dean  and  Dalrymple,  2002).  The  equilibrium  beach 
profile  relates  the  long-term  shape  of  the  profile  elevation  of  the 
beach  profile,  y ,  to  distance  offshore,  x, 

x  =  Ay2'3.  (9) 

During  an  erosion  event,  the  eroded  sand  is  transported  offshore 
and  deposited,  thereby  decreasing  depths  offshore  such  that  sediment 
volume  is  conserved.  Avalanching  of  the  profile  is  initiated  if  the  slope 
is  greater  or  equal  to  a  user-specified  avalanching  angle  with  a  default 
value  of  30°. 

3.3.  Runup  overwash  and  inundation  overwash 

Overwash  is  any  wave  uprush  which  passes  over  the  “crown,”  or 
crest  of  the  barrier  beach  (Leatherman,  1979,  p.  3).  Of  relevance  to  this 
study  is  the  magnitude  of  the  morphologic  feature  created  by 


overwash  and  deposited  on  the  bayside  of  the  crest,  called  “washover,” 
“washover  deposit”  or  “washover  fan”  (Leatherman,  1979,  p.  2).  The 
frequency  and  magnitude  of  overwash  depend  on  long-term  condi¬ 
tions,  such  as  storm  climatology,  relative  sea  level  rise,  and  sediment 
supply.  Overwash  and  the  resulting  washover  deposit  are  one  of  the 
mechanisms  through  which  the  barrier  island  migrates  towards  the 
bay  (across  shore).  Two  modes  of  overwash  are  simulated  in  the 
model:  runup  overwash  and  inundation  overwash.  Runup  overwash 
occurs  if  the  island  is  not  submerged,  and  washover  is  caused  by  the 
uprushing  wave  bore.  Inundation  overwash  occurs  when  the  storm 
surge  level  and  wave  setup  exceed  the  elevation  of  the  barrier  island 
crest,  and  the  barrier  island  is  submerged  (Donnelly  et  al.,  2009). 

The  overwash  transport  rate  over  the  beach  crest  due  to  runup 
overwash  per  unit  length  of  beach,  qDR(t),  can  be  described  as 
(Donnelly  et  al.,  2009), 

|  (10) 

for  0  <  zR{t)  and  S(t)  +  d(t)  <  bh(t) 

where  I<R  is  a  calibration  coefficient  that  accounts  for  sediment 
concentration  and  properties  of  the  wave  bore,  zR{t)  is  the  elevation  of 
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the  runup,  R(t),  relative  to  the  dune  crest  elevation,  bh(t),  and  d(t)  is 
the  local  water  depth  (Fig.  3b).  For  calculations  herein,  I<R  was  set 
equal  to  0.005  as  recommended  by  Donnelly  et  al.  (2009). 

The  two-percent  runup,  Ru2%(t),  is  calculated  as  (Hughes,  2004), 


^u2  %(f) 


4.4(S(t)  +  d(t))  tan  fJ(t)0,70 


MF(t) 

pg(S(t)  +  d(t))2 


1/2 


for  —  <  tan  (5(t)  < 


1 

5 


(11) 


in  which  p  is  the  density  of  water,  and  the  maximum  dimensionless 
depth-integrated  wave  momentum  flux  per  unit  width,  MF(t),  is 


MF(t) 

=  A>(t) 

'm  +  d(t»' 

p sm  +  d(t))\ 

max 

g7p(t)2 

where  A0(t)  =  0.6392 


(  Hmoffl  V°256 
V(S(t)  +  d(t))J 


(12) 


and  A,(t)  =  0.1804 


(  Hmo(t)  \ 


-0.391 


The  zeroth-moment  deep-water  wave  height  is  Hmo(t),  with 
associated  peak  period,  Tp(t).  If  wave  height  is  randomly  generated 
about  a  user-specified  mean,  2D  MCO  calculates  the  peak  wave  period, 
Tp(t),  associated  with  Hmo(t)  using  a  relationship  from  Bretschneider 
(1966)  based  on  hurricanes  which  has  been  modified  for  peak  wave 
period, 

yt)  (sec)  =  2.13v/Hmo(t)(ft)  =  3.86  vTWtXm).  (13) 


The  transport  rate  over  the  beach  crest  per  unit  width  of  beach  due 
to  inundation  overwash,  qDI(t),  is  given  by  (Donnelly  et  al.,  2009)  as, 

qDI(t)  =  2K|y/2gzR(t)3/2  +  qD  R(t), 

(14) 

0  <zR{t)  and  S(t)  +  d(t)>bh(t) 

in  which  I<{  is  an  empirical  coefficient,  and  zR(t)  is  as  defined  pre¬ 
viously  (Fig.  3c).  For  calculations  herein,  Kx  was  set  to  0.005  (Donnelly 
et  al.,  2009). 

Transport  in  the  swash  zone  is  calculated  as  (Larson  et  al.,  2004), 
<?sw(t)  =  ^sw2\/2gRU2%(t)3/2(tan(isw-tan(ieq)  (15) 


where  I<sw  is  an  empirical  coefficient,  set  equal  to  0.0016  for 
calculations  herein,  /3SW  is  the  local  slope  in  the  swash  zone,  and 
jBeq  is  the  equilibrium  slope  calculated  using  equilibrium  profile 
concepts, 


“’eq 


-a(; 


Y2  /  _ y‘ 

Adune  A< 


,2/3 

sw 


y  ^dune 


-x* 


(16) 


Location  of  the  dune  and  swash  are  given  by  xdUne  and  xsw, 
respectively. 


3.4.  Consolidation 


The  original  relationship  for  calculating  soil  consolidation  was 
derived  by  Terzaghi  (1925,  1943)  and  Terzaghi  and  Peck  (1967). 
Primary  consolidation  is  defined  as  the  process  during  which  excess 
pore  water  pressure  is  dissipated  from  a  soil  particle  matrix,  and  is 


derived  from  fundamental  hydraulic  principles.  The  assumptions  for 
one-dimensional  consolidation  theory  are:  (1)  a  fully-saturated 
sediment  system;  (2)  unidirectional  flow  of  water;  (3)  one-dimen¬ 
sional  compaction  occurring  in  the  opposite  direction  of  flow;  (4)  a 
linear  relationship  between  the  change  in  sediment  volume  and  the 
applied  pressure  (linear  small-strain  theory);  and  (5)  validity  of 
Darcy's  Law,  which  states  that  the  specific  discharge  (flow  rate  per 
area)  through  a  porous  medium  is  equal  to  the  hydraulic  gradient 
times  the  hydraulic  conductivity  (Yong  and  Warkentin,  1966; 
Hornberger  et  al.,  1998).  For  one-dimensional  vertical  flow,  the 
maximum  consolidation,  zc,  for  a  given  loading  p  greater  than  an 
initial  loading,  px,  can  be  calculated  as, 


where  z0  is  the  thickness  of  compressible  sediment;  Cx  is  the  value  of 
the  compression  index  that  is  determined  experimentally  from  a 
Casagrande  consolidation  test  (Casagrande  and  Fadum,  1940;  Horn¬ 
berger  et  al.,  1998);  and  e0  is  the  initial  void  ratio,  equal  to  the  volume 
of  voids  divided  by  the  volume  of  solids  as  averaged  over  the  thickness 
z0.  The  parameter  z0  can  be  estimated  from  sediment  core  data  or 
depositional  maps  that  represent  the  thickness  of  soft  sediment.  For 
the  Mississippi  River  system,  Kulp  et  al.  (2002)  analyzed  data  from 
more  than  800  boreholes  and  mapped  the  topstratum  lithosome, 
which  represents  fine-grained  deposition  in  fluvial,  deltaic,  and  shelf 
environments  that  overlies  coarser-grained  substratum.  Kulp  et  al. 
determined  that  the  maximum  thickness  of  topstratum  sediment  is 
approximately  120  m  in  vicinity  of  the  modern  Balize  depocenter,  and 
ranges  from  approximately  10  to  30  m  in  vicinity  of  the  modern 
barrier  islands  (Kulp  et  al.,  2002,  their  Fig.  7).  In  the  applications 
discussed  herein,  z0  did  not  vary  spatially,  but  non-homogeneous 
variation  in  the  subsurface  characteristics  could  be  incorporated  in 
site-specific  application  of  the  model. 

The  magnitudes  of  Cx  and  px  vary  depending  on  whether  p  is  greater 
or  less  than  the  pre-consolidation  loading,  pc,  which  can  be  determined 
from  a  Casagrande  consolidation  test.  The  magnitude  of  the  pre¬ 
consolidation  stress  is  important  because  it  separates  soils  that  are  over¬ 
consolidated  (i.e.,  these  soils  have  experienced  a  greater  load  at  some 
time  in  their  past)  from  those  that  are  under-consolidated  (i.e.,  the 
present  loading  is  the  maximum  that  has  occurred).  Loading  greater 
than  the  pre-consolidation  stress  will  result  in  greater  rates  of 
consolidation  than  have  previously  occurred.  Fig.  4  shows  results  of 
a  Casagrande  consolidation  test  conducted  for  a  sediment  sample  at 
12.5-13.1  m  depth  from  Chaland  Headland,  Louisiana. 

For  the  example  shown  in  Fig.  4,  if  the  loading  p  is  greater  than 
the  pre-consolidation  stress,  pc,  then  px  =  pc  =  7900  kg/m2,  and  Cx  = 
Cc  =  0.4  in  Eq.  (17).  If  p  is  less  than  the  pre-consolidation  stress, 
then  px  is  equal  to  the  initial  stress,  pi  =  660  kg/m2,  and  Cx  =  Cd  = 
0.125  in  Eq.  (17). 

Terzaghi's  ( 1 943 )  time-dependent  relationship  for  consolidation  is 


du  _  d2u 
dt  ~  CmW 


(18) 


where  u  is  pore  water  pressure  in  excess  of  hydrostatic  pressure,  t  is 
elapsed  time,  cvx  is  a  property  of  the  compressible  sediment,  called  the 
coefficient  of  consolidation  and  also  determined  during  consolidation 
testing,  and  z  is  the  vertical  coordinate  with  the  origin  at  the  initial 
sediment  surface  (Fig.  5).  If  the  loading  is  greater  than  the  pre¬ 
consolidation  stress,  then  cvx  =  cvc ;  otherwise  cvx  =  cvi  (Fig.  4). 

The  proportion  of  the  initial  pore  water  pressure  remaining  at  any 
time,  M(t),  can  be  expressed  as, 


M(t) 


1  Po  iL  dz  = 
z0J0  u0  e0— ef 


(19) 
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Applied  Axial  Stress,  p„  kg/m2 


in  which  u0  is  the  initial  pore  water  pressure,  e(t)  is  the  average  void 
ratio  at  any  time,  and  ef  is  the  final  average  void  ratio  corresponding  to 
the  consolidation  test  results  for  the  portion  of  the  curve  less  than  or 
greater  than  the  pre-consolidation  stress.  The  variable  M(t)  varies 
between  1  and  0,  at  time  t=  0  and  infinity,  respectively.  The 
proportion  of  vertical  consolidation  that  occurs  at  any  time  can  also 
be  expressed  as, 

«>  =  ««  m 

Combining  Eqs.  (19)  and  (20)  gives, 

z(x,t)  =  zc(x,t)(l-M(t))  (21) 


where  M(t)  can  be  expressed  as  (Dean,  2002,  p.  119) 


M(t) 


00  l(2n  l)n]2cvxt / 4zq 

n?i  (2n— l)2n2 


(22) 


Fig.  6  shows  a  flowchart  for  the  consolidation  routine  in  the  2D 
MCO  model. 


The  consolidation  routine  was  evaluated  with  information  pre¬ 
sented  by  Blum  et  al.  (2008)  in  which  they  discussed  uplift  and 
subsidence  of  the  Mississippi  River  Delta  over  the  past  30,000  years. 
Although  the  spatial  and  temporal  scales  applied  by  Blum  et  al.  for  the 
Mississippi  River  Delta  are  much  greater  than  for  a  barrier  island,  the 
magnitude  of  consolidation  as  a  function  of  the  thickness  and 
associated  weight  of  sediment  are  commensurate.  Specific  to  this 
study,  Blum  et  al.  discuss  the  gradual  deposition  of  approximately 
40  m  of  deltaic  sediment  in  the  vicinity  of  New  Orleans  (90.5° 
longitude)  from  11,500  to  4,000  years  ago,  and  no  additional 
deposition  from  4000  years  ago  to  the  present  as  the  river's  depo- 
center  had  moved  further  downstream.  Blum  et  al.'s  calculations  with 
a  three-dimensional  visco-elastic  consolidation  model  indicated  that 
approximately  5.9  m  surface  deflection  (the  net  of  consolidation  plus 
uplift)  occurred  over  the  past  10,000  years  (see  their  Fig.  3b  and  c, 
reproduced  below  in  Fig.  7b  and  c).  As  the  river  incised  the  Lower 
Mississippi  River  Valley  with  meltwater,  the  removal  of  sediment 
created  an  uplift  of  the  surface  from  30,000  to  9500  years  ago.  This 
loading  cycle  (linear  deposition  of  40  m  thickness  from  11,500  to 
4000  years,  then  no  deposition  from  4000  years  to  present)  was 
programmed  into  the  consolidation  code  with  the  same  sediment- 
water  density  as  discussed  by  Blum  et  al.  (1800  kg/m3)  and  indicated 
a  value  of  6.7  m  consolidation  over  the  past  11,500  years  (Fig.  7).  The 
value  calculated  by  the  consolidation  routine  developed  herein  is  15% 
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Fig.  5.  Definition  sketch  for  consolidation  relationship. 


8 


J.D.  Rosati  et  al.  /  Marine  Geology  271  (2010)  1-16 


.! 

1 

OJ 

LD 

PS 


Initial  subsurface 


Consolidated  subsurface 


Initial  profile,  y(x,1)' 

Time-dependent 

consolidation 

ZtotM 


Profile  at  later  time,  y(x,t) 


V  Mean  Sea 


=  Leve 


x,  Alongshore  distance 


Calculation  ax 
cells 


Fig.  6.  Flowchart  for  2D  MCO  consolidation  calculations. 


greater  than  Blum  et  al.'s  results.  However,  Blum  et  al.'s  value  includes 
approximately  0.5  m  that  occurred  at  the  end  of  the  uplift  period. 
Adjusting  the  consolidation  calculations  for  the  uplift  results  in  an 
estimate  of  6.2  m,  within  5  %  of  Blum  et  al.'s  results.  This  comparison  is 
considered  a  validation  of  the  consolidation  routine. 

4.  Model  comparison  with  field  data 

To  fully  calibrate  and  validate  2D  MCO,  an  initial  beach  profile  and 
subsurface  elevation,  similar  data  several  years  to  decades  later, 
sediment  characteristics  of  the  overlying  sand  beach  and  compress¬ 
ible  subsurface,  and  the  forcing  hydrodynamics  during  that  time 
period  would  be  required.  A  complete  search  of  existing  databases 
and  literature  was  conducted  and  identified  only  a  data  set  from 


Virginia  as  sufficient  for  comparison  with  model  applications.  The 
term  “comparison”  is  introduced  here  to  describe  application  of  the 
model  to  measurements  in  which  several  unknowns  were  estimated. 
In  this  comparison,  the  model  input  conditions  were  likely  only  one  of 
many  sets  that  could  be  developed  and  executed  in  the  model  to 
adequately  reproduce  the  measured  response.  Although  the  compar¬ 
ison  to  Virginia  data  is  not  an  ideal  calibration  and  validation  of  the 
model,  it  lends  credibility  to  2D  MCO  applications. 

4.1.  Virginia  data 

In  a  study  of  Virginia  barrier  islands,  Gayes  (1983)  surveyed  the 
barrier  and  beach  profile  and  collected  sediment  cores  across  three 
migrating  barrier  island  systems  that  overlie  a  compressible  peat  and 
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Fig.  7.  Comparison  of  consolidation  routine  (shown  in  Fig.  7a)  to  calculations  from  Blum  et  al.  (2008)  (reproduced  in  7b  and  7c,  their  Figs.  3b  and  3c). 


bay  sediment  substrate:  Assawoman  Island,  Metomkin  Island,  and 
Wallops  Island,  Virginia,  U.S.A  (Figs.  8-11).  Wallops  Island  represents 
a  site  disturbed  by  construction  and  adjacent  inlet  processes  and  is  not 
discussed  further  (see  Gayes,  1983  for  more  detail).  When  adjusted 
for  eustatic  sea  level  rise  of  2  mm/yr  (Douglas,  1992;  Peltier,  1998) 
over  a  40-year  period,  these  barrier  island  systems  experienced 
consolidation  of  the  substrate  between  0.9  and  1.1  m  as  they  migrated 
landward  at  rates  ranging  from  3.8  to  4.8  m/yr.  The  elevations  of  these 
islands  were  approximately  1.2  to  2.6  m  relative  to  Mean  High  Water 
(MHW),  with  the  maximum  thickness  of  sand  overlying  the  substrate 
ranging  from  1.4  to  3.5  m.  If  consolidation  of  the  underlying  substrate 


had  not  occurred,  the  sand-clay/silt  interface  would  lie  at  approxi¬ 
mately  the  zero  MHW  line.  Void  ratios  of  the  back-barrier  sediment 
are  greater  than  those  of  the  clay  and  silt  underlying  the  sandy  barrier 
island,  reinforcing  the  interpretation  of  consolidation  due  to  the 
loading  of  the  island.  Measurements  are  summarized  in  Table  1. 

The  goals  of  model  comparisons  with  the  Virginia  data  were  to  test 
whether  2D  MCO  could  reproduce  the  cross-shore  shape  and  elevation  of 
the  profiles  and  magnitude  and  shape  of  the  consolidated  surface. 
Unknowns  in  this  process  were  the  initial  profile  shape,  storm  wave  and 
surge  conditions  during  the  simulation  period,  the  magnitude  ofz0  in  the 
consolidation  calculations  ( Eqs.  ( 1 7 ),  ( 1 9 ),  and  ( 22 ) ),  and  results  from  a 


Assawoman  Island  Transect  1 
Assawoman  Island  Transect  2 


Atlantic  Ocean 


Metomkin  Island  Transect 


Fig.  8.  Location  of  Virginia  data  (adapted  from  Gayes  1983). 
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Fig.  9.  Assawoman  Island,  Virginia,  cross-section  1  (adapted  from  Gayes  1983). 


Casagrandeconsolidationtestofthesubsurfaceorback-barriersediment 
cores  ( eG,  Cx,  and  px  in  Eq.  ( 1 7 ) ) .  Values  of  cvx  ( Eq.  ( 22 ) )  were  available  for 
Metomkin  Island  at  two  locations,  and  the  average  value  of  a  back-barrier 
sample  equal  to  2.5  m2/yr  was  applied  for  the  Metomkin  and  Assawoman 
sites.  In  the  absence  of  Casagrande  consolidation  test  data,  values  for 
ChalandHeadland,Louisiana(Fig.4)wereadoptedfortheVirginiabarrier 
islands. 


All  islands  have  marsh  on  the  back-barrier,  and  the  2D  MCO  model 
was  modified  to  allow  the  back-barrier  marsh  surface  to  maintain 
elevation  with  MHW  during  the  simulation  period.  The  mean,  spring, 
and  neap  tidal  amplitudes  for  Metomkin  are  0.55,  0.65,  and  0.43  m, 
respectively  (Byrnes  and  Gingerich,  1987).  Based  on  Gayes'  descrip¬ 
tion  that  overwash  was  “extensive”  for  the  Metomkin  site  (at  1.2  m 
MHW  dune  elevation)  and  “infrequent”  for  the  Assawoman  sites  (at 
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Fig.  10.  Assawoman  Island,  Virginia,  cross-section  2  (adapted  from  Gayes  1983). 
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Fig.  11.  Metompkin  Island,  Virginia  (adapted  from  Gayes  1983). 


1.4  and  2.6  m  MHW  dune  elevations),  it  is  speculated  that  a  dominant 
storm  water  level  elevation  (including  surge,  tide,  wave  setup  and 
runup)  has  been  between  1.2  and  2.6  m  MHW  for  this  region  of  the 
coast  during  the  period  of  measurement.  Wave  Information  Study 
(WIS)  data  for  Station  179,  located  in  18-m  water  depth  at  latitude 
37.75°  and  longitude  75.33°,  indicated  that  the  mean  Hmo  =  0.9± 
0.6  m,  with  associated  peak  period  Ip  =  6±2.8  s  (WIS,  2008).  From 


these  discussions  and  in  the  absence  of  other  data,  it  is  assumed  that  a 
reasonable  value  for  total  storm  water  elevation  at  the  site  ranges 
from  1.2  to  2.6  m  MHW,  with  maximum  deep-water  significant  wave 
height  approximately  1.5  m. 

For  comparison  with  2D  MCO,  variations  in  storm  surge,  offshore 
wave  height,  and  tidal  elevation  (randomly  generated  within  the 
mean  tide  amplitude  =+/— 0.55  m)  were  applied  to  force  the  model 


Table  1 

Data  from  Virginia. 
Adapted  from  Gayes,  1983. 


Site1 

Dune  elev 
(m  MHW) 

Max  sand  below 
(m  MHW) 

Sand  vol 

(m3/m) 

Mig  rate2 
(m/yr) 

Relative 

back-barrier  sed 

Relative  freq 
of  overwash 

Relative 
sand  supply 

Max  consol  due 
to  weight  of  sand3 

(m) 

Al  (Fig.  8) 

1.4 

2.2 

225 

3.75 

Finest;  organic-rich 

Infrequent 

Greatest 

1.1 

A2  (Fig.  9) 

2.6 

3.5 

343 

4.05 

Coarsest;  greatest  energy 

Infrequent 

Intermediate 

1.0 

M  (Fig.  10) 

1.2 

1.4 

194 

4.79 

Medium 

Extensive 

Lowest 

1.2 

W= Assawoman  1,  A2=Assawoman  2,  M=Metomkin. 

2Over  a  40-year  period. 

3Adjusted  to  remove  sequestering  of  sand  due  to  eustatic  sea  level  rise  =  2  mm/yr  (revised  from  original  calculation  by  Gayes  (1983)  in  which  he  applied  1  mm/yr). 


Table  2 

Summary  of  2D  MCO  applications  to  Virginia  Barrier  Islands  (100  year  calculation). 


2D  MCO  input  parameters 

Calculations  and  measurements 

Initial  dune  elevation 
(m  MHW) 

Island  width 

(m) 

Hmo 

(m) 

Avg  total  storm  water  level1 
(m  MHW) 

Zo 

(m) 

Total  consolidated  volume  Dune  elevation 

(m3/m  (%V))  (m  MHW) 

Maximum 
consolidation 
(m  MHW) 

Calc  Meas  Calc  Meas 


Assawoman  cross-section  1 

5.6  400  1.5  1.5 

Assawoman  cross-section  2 

7.6  400  1.5  2.0 

Metomkin 

7.6  400  1.5  1.8 


10 

460  (50%) 

1.4 

1.4 

1.1 

1.1 

9 

577 (39%) 

2.4 

2.6 

1.2 

1.0 

11 

739  (56%) 

1.3 

1.2 

1.3 

1.2 

Average  total  storm  water  level  includes  storm  surge,  wave  runup,  tide,  and  relative  sea  level  rise  over  100-year  simulation  period. 
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Fig.  12.  Summary  of  2D  MCO  calculations  for  Assawoman  cross-section  1  simulation  as  measured  at  the  dune  crest. 


over  a  100-year  simulation  period,  which  was  a  sufficient  duration  to 
allow  the  barrier  island  and  subsurface  to  equilibrate  to  the  forcing 
and  loading  conditions.  In  the  model  simulations,  eustatic  sea  level 
rise  was  2  mm/yr  (Douglas,  1992;  Peltier,  1998)  and  water  level 
increased  each  year  with  sea  level  rise. 

An  iterative  process  was  applied  in  the  modeling  exercise, 
modifying  the  starting  dune  crest  elevation,  barrier  island  width, 
wave  and  surge  conditions,  and  thickness  of  compressible  sediment  z0 
such  that  the  measured  profile  cross-sections  and  cross-shore 
magnitude  of  consolidation  were  approximated  after  the  simulation 
period.  Eight  thousand  1-m  cross-shore  calculation  cells,  and  up  to 
thirty  1-m  vertical  cells  were  simulated  with  the  Matlab©  code,  which 
took  approximately  2  to  5  min  for  100-year  simulation  on  a  1.86  GHz 


personal  computer.  Input  variables  and  a  comparison  between 
measured  and  calculated  output  for  simulations  that  agreed  best 
with  the  known  information  about  the  sites  are  summarized  in 
Table  2.  Fig.  12  summarizes  results  from  the  model  simulations  as 
measured  at  the  dune.  Simulation  profiles  calculated  after  100  years, 
as  well  as  the  subsurface  at  25,  50,  and  100  years  are  shown  for  each 
site  in  Figs.  13-15. 

Model  calculations  approximate  the  general  shape  of  the  beach  and 
subsurface  elevations,  with  the  maximum  dune  crest  elevation  and 
minimum  consolidated  subsurface  elevation  agreeing  within  0.2-m  of 
measured  values.  The  calculated  profile  poorly  represents  the  measured 
beach  profile  seaward  of  the  dune,  likely  because  erosion  and  overwash 
processes  are  only  represented  during  storm  events  when  water  levels 


2050  2100  2150  2200  2250  2300  2350  2400 


Cross-Shore  Distance,  m 

Fig.  13.  Comparison  of  2D  MCO  calculations  with  measured  beach  and  subsurface  profiles  for  Assawoman  Island  cross-section  1. 
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Fig.  14.  Comparison  of  2D  MCO  calculations  with  measured  beach  and  subsurface  profiles  for  Assawoman  Island  cross-section  2. 


are  near  the  dune  crest  elevation.  The  shape  of  the  calculated  subsurface 
agrees  with  the  general  trend  of  the  measured  subsurface  at  the  25  and 
50-year  durations.  The  shape  of  the  calculated  100-year  subsurface 
indicates  additional  compression  from  washover  sand  on  the  back- 
barrier  marsh,  which  is  not  evident  in  the  measured  subsurface  but 
becomes  more  likely  over  longer  time  periods. 

5.  Implications  of  compressible  substrate 

In  this  section,  three  sets  of  2D  MCO  simulations  are  presented  to 
compare  evolution  of  barrier  islands  overlying  a  compressible 
substrate  (BIC)  to  those  overlying  a  stable  substrate  (BIS).  The  first 
series  evaluated  how  compressibility  of  the  substrate  modified  cross¬ 
shore  barrier  island  migration  and  morphology  change  over  a  50-year 
simulation  period.  The  analysis  varied  the  thickness  of  compressible 
sediment  (z0)  for  a  barrier  island  with  initial  dune  elevation  of  3.5  m 
relative  to  mean  sea  level  (MSL),  2200  m  width  at  MSL,  1-m  average 
storm  wave  height,  and  -him  average  surge  relative  to  MSL.  Five 


magnitudes  of  z0  were  evaluated:  0  (stable  substrate),  5,  10,  20,  and 
25  m,  with  consolidation  parameters  cvi  =  cvc  =  2.5  m2/yr.  Simulations 
terminated  if  the  crest  of  the  barrier  island  became  submerged  below 
MSL.  Results  are  shown  in  Fig.  16. 

Larger  values  of  z0  increased  the  magnitude  of  maximum  consoli¬ 
dation,  the  rate  of  dune  lowering,  and  decreased  the  duration  of  the 
simulation.  The  migration  rate  increased  as  substrate  thickness 
increased  from  0  to  10  m  then  was  constant  for  greater  substrate 
thicknesses.  The  constant  migration  for  substrate  thickness  greater  than 
20  m  reflects  the  tendency  for  BIC  to  drown  in  place  rather  than  migrate, 
as  the  consolidation  process  dominates  the  island  evolution.  This 
process  of  island  lowering,  followed  by  breakup  and  submergence  has 
been  observed  for  barrier  islands  with  a  limited  sand  supply  that  overlie 
thick  deltaic  deposits  in  Louisiana  (Penland  and  Boyd,  1981;  Kuecher, 
1994;  Roberts  et  al.,  1994;  McBride  et  al.,  1995;  Kulp  et  al.,  2002). 

The  next  series  of  simulations  evaluated  whether  a  sufficient 
supply  of  littoral  sand  would  increase  the  migration  rates  for  BIC.  A 
longshore  sand  transport  (LST)  sub-module  was  coupled  with  2D 


Fig.  15.  Comparison  of  2D  MCO  calculations  with  measured  beach  and  subsurface  profiles  for  Metomkin  Island. 
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Fig.  16.  Influence  of  thickness  of  compressible  substrate  on  sandy  barrier  island  evolution. 


MCO  to  provide  a  net  source  of  sand  to  the  island  (arbitrarily  set  to  a 
maximum  of  75  m3/m/yr)  if  it  was  submerged  below  MSL  during  the 
previous  storm  season.  The  source  of  LST  implies  that  the  rate  of  sand 
transport  provided  to  the  island  is  greater  than  that  removed  from  the 
island;  if  LST  increases  with  distance  along  the  island  or  diverges,  then 
LST  represents  a  net  sink  (or  loss)  to  the  island.  In  the  simulations,  if 
the  source  of  LST  was  insufficient  to  maintain  the  island  above  MSL, 
the  simulation  terminated.  Simulations  were  set  for  a  99-year 
duration  with  the  same  island  and  storm  parameters  as  in  the 
previous  simulations,  and  the  thickness  of  compressible  substrate 
zo  =  0,  10,  15,  and  20  m  with  consolidation  parameters  cvi  =  cvc  = 
2.5  m2/yr.  Results  of  these  simulations  are  shown  in  Fig.  17. 

Net  LST  rates  required  to  maintain  the  islands  through  the 
migration  process  were  equal  to  4m3/m/yr  for  BIS  and  ranged 
between  31  and  73  m3/m/yr  for  BIC.  Migration  rates  of  BIC  were  3  to 
11  times  greater  than  for  BIS,  if  a  sufficient  source  of  net  LST  was 
available  to  replenish  the  island.  Nevertheless,  the  LST  rate  was  not 
sufficient  to  maintain  BIC  with  z0  =  15  and  20  m  for  the  full  99  years; 
simulations  terminated  after  60  and  43  years,  respectively.  Thus,  for 
this  set  of  simulation  parameters,  deltaic  deposits  greater  than  15  m 


thickness  continue  to  dominate  morphologic  change  and  cause  island 
submergence  even  with  a  net  source  of  LST. 

The  final  series  evaluated  whether,  in  the  absence  of  a  sufficient 
LST  sand  source,  the  lifetime  of  BICs  with  varying  dune  crest 
elevations  was  reduced  as  compared  to  BIS.  Simulations  varied  initial 
dune  crest  elevation  from  2  to  4.5  m  MSL  for  a  50-year  simulation 
period,  and  are  summarized  in  Table  3.  All  BIC  durations  were  less 
than  those  of  BIS  except  for  the  largest  barrier  island  (dune  crest 
elevation  =  4.5  m  MSL).  For  these  cases,  the  island  provided  a 
sufficient  elevation  to  reduce  washover  as  evidenced  by  similar 
migration  rates  for  both  BIC  and  BIS.  The  simulations  suggest  that  the 
lifespan  of  BIC  are  reduced  as  compared  to  BIS,  except  for  BIC  that  are 
of  sufficient  elevation  to  reduce  washover  over  the  time  periods  of 
consideration  and  for  prevailing  storm  conditions  at  the  site. 

6.  Conclusions 

This  study  developed  a  cross-shore  model  of  sandy  barrier  island 
erosion  and  overwash,  including  time-dependent  consolidation  of  the 
underlying  compressible  substrate  as  a  function  of  the  thickness  and 


Fig.  17.  Evolution  of  islands  overlying  a  compressible  substrate  with  a  net  source  of  sand  to  replenish  sand  budget. 
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Table  3 

Lifespan  and  migration  rates  for  BIC  and  BIS  for  a  50-year  simulation. 


Initial  dune  elevation 
(m  MSL) 

Duration 

(yr) 

Migration 

(m/yr) 

BIC 

BIS 

BIC 

BIS 

2 

15 

50 

205 

10.8 

2.5 

34 

50 

84 

7.4 

3.5 

46 

50 

78 

8.6 

4.5 

50 

50 

7.3 

8.2 

sediment  characteristics  of  the  underlying  substrate,  and  weight  of 
the  barrier  island.  The  model  is  intended  for  applications  spanning 
years  to  decades  and  was  applied  herein  for  simulations  of  100  years 
in  duration  to  test  the  hypothesis  that  substrate  characteristics  modify 
the  cross-shore  evolution  (migration  and  volumetric  losses)  of  a 
barrier  island.  The  hypothesis  was  formulated  to  evaluate  whether  the 
consolidation  process  accelerates  long-term  morphologic  response  for 
sandy  barrier  islands  overlying  poorly-consolidated  sediment,  such  as 
deltaic,  bay,  estuarine,  and  peat  deposits.  Data  to  validate  the  model 
were  limited  to  profile  and  subsurface  elevations  for  three  cross-shore 
transects  of  two  barrier  islands  in  Virginia.  Reasonable  values  for  other 
input  parameters  were  estimated  to  compare  with  model  calcula¬ 
tions.  After  calibration,  the  model  was  able  to  reproduce  the 
maximum  dune  crest  and  minimum  subsurface  elevation  within 
0.2-m  of  the  measured  values  in  a  100-year  simulation. 

Application  of  the  model  validated  the  hypothesis  through 
simulations  of  cross-shore  barrier  island  morphologic  change.  For  a 
compressible  substrate  thickness  from  10  to  30  m,  such  as  beneath  the 
modern  Louisiana  barrier  islands  (Kulp  et  al.,  2002),  model  simulations 
indicated  that  substrate  thickness  was  non-linearly  related  to  the 
consolidation  rate  and  loss  of  dune  crest  elevation.  For  substrate 
thickness  greater  than  10  m,  model  calculations  indicated  a  tendency 
for  barrier  islands  to  drown  in  place  as  the  process  of  subsurface 
consolidation  dominated  the  morphologic  response.  With  a  sufficient 
source  of  sand  to  the  island,  migration  rates  for  islands  overlying  a 
compressible  substrate  were  3  to  1 1  times  greater  than  for  the  same 
island  over  a  stable  substrate.  Without  a  source  of  sand,  islands  over- 
lying  a  compressible  substrate  became  submerged  more  rapidly  than 
within  a  stable  substrate.  The  exception  was  for  a  barrier  island  with 
sufficient  elevation  such  that  overwash  and  inundation  did  not  occur. 

Barrier  islands  overlying  a  compressible  substrate  incur  an 
additional  volumetric  loss  due  to  the  consolidation  process  as  a 
function  of  the  magnitude  of  loading  applied  to  the  substrate  and 
duration  of  the  loading.  Because  of  consolidation  of  the  substrate, 
islands  that  migrate  and  overwash  experience  additional  losses  due  to 
consolidation,  as  compared  to  islands  that  are  stable.  Barrier  islands 
overlying  a  consolidating  substrate  are  more  likely  to  have  ( 1 )  reduced 
dune  elevations  because  of  consolidation,  (2)  overall  volumetric  loss 
to  fill  in  compressed  regions  outside  the  footprint  of  the  island,  and 
(3)  increased  overwash  and  migration  after  the  dune  reaches  a  critical 
elevation  with  respect  to  the  total  water  elevation  of  the  prevalent 
storm  conditions.  In  effect,  the  consolidation  process  decreases  the 
return  period  (increases  the  frequency)  of  the  prevailing  storm  con¬ 
ditions.  Numerical  calculations  with  the  model  illustrated  how  consol¬ 
idation  modifies  profile  response  through  lowering  of  the  dune 
elevation  and  increasing  the  potential  for  overwash  and  migration. 
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